1 

Distributed Transmit Beamforming using 

Feedback Control 

R. Mudumbai, Student Member, IEEE, J. Hespanha, Senior Member, IEEE, U. Madhow, Fellow, IEEE 

and G. Barriac, Member, IEEE 



Abstract — A simple feedback control algorithm is pre- 
sented for distributed beamforming in a wireless network. 
A network of wireless sensors that seek to cooperatively 
transmit a common message signal to a Base Station (BS) 
is considered. In this case, it is well-known that substan- 
tial energy efficiencies are possible by using distributed 
beamforming. The feedback algorithm is shown to achieve 
the carrier phase coherence required for beamforming 
in a scalable and distributed manner. In the proposed 
algorithm, each sensor independently makes a random 
adjustment to its carrier phase. Assuming that the BS is 
able to broadcast one bit of feedback each timeslot about 
the change in received signal to noise ratio (SNR), the 
sensors are able to keep the favorable phase adjustments 
and discard the unfavorable ones, asymptotically achiev- 
ing perfect phase coherence. A novel analytical model is 
derived that accurately predicts the convergence rate. The 
analytical model is used to optimize the algorithm for 
fast convergence and to establish the scalability of the 
algorithm. 
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I. Introduction 

Energy efficient communication is important for wire- 
less ad-hoc and sensor networks. We consider the prob- 
lem of cooperative communication in a sensor network, 
where there are multiple transmitters (e.g., sensor nodes) 
seeking to transmit a common message signal to a distant 
Base Station receiver (BS). In particular, we investigate 
distributed beamforming, where multiple transmitters 
coordinate their transmissions to combine coherently 
at the intended receiver. With beamforming, the sig- 
nals transmitted from each antenna undergo constructive 
interference at the receiver, the multiple transmitters 
acting as a virtual antenna array. Thus, the received 
signal magnitude increases in proportion to number of 
transmitters N, and the SNR increases proportional to 
TV 2 , whereas the total transmit power only increases 
proportional to N. This iV-fold increase in energy ef- 
ficiency, however, requires precise control of the carrier 
phases at each transmitter in order that the transmitted 
signals arrive in phase at the receiver. In this paper, 
we propose a feedback control protocol for achieving 
such phase coherence. The protocol is based on a fully 
distributed iterative algorithm, in which each transmitter 
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independently adjusts its phase by a small amount in a 
manner depending on a single bit of feedback from the 
BS. The algorithm is scalable, in that convergence to 
phase coherence occurs in a time that is linear in the 
number of cooperating transmitters. 

Prior work on cooperative communication mainly fo- 
cuses on exploiting spatial diversity for several wireless 
relaying and networking problems [1], [2]. Such dis- 
tributed diversity methods require different transmitters 
to transmit information on orthogonal channels, which 
are then combined at the receiver. The resulting diver- 
sity gains could be substantial in terms of smoothing 
out statistical fluctuations in received power due to 
fading and shadowing environments. However, unlike 
distributed beamforming, distributed diversity does not 
provide a gain in energy efficiency in terms of average 
received power, which simply scales with the transmitted 
power. On the other hand, the coherent combining of 
signals at the receiver due to distributed beamforming 
also provides diversity gains. 

Recent papers discussing potential gains from dis- 
tributed beamforming include [3], which investigates the 
use of beamforming for relay under ideal coherence 
at the receiver, and [4], which shows that even partial 
phase synchronization leads to significant increase in 
power efficiency in wireless ad hoc networks. The beam 
patterns resulting from distributed beamforming using 
randomly placed nodes are investigated in [5]. However, 
the technical feasibility of distributed beamforming is 
not investigated in the preceding papers. In our prior 
work [6], [7], we recognized that the key technical 
bottleneck in distributed beamforming is carrier phase 
synchronization across cooperating nodes. We presented 
a protocol in which the nodes first establish a common 
carrier phase reference using a master-slave architecture, 
thus providing a direct emulation of a centralized multi- 



antenna system. This is a challenging problem, because 
even small timing errors lead to large phase errors at the 
carrier frequencies of interest. Once phase synchroniza- 
tion is achieved, reciprocity was proposed as a means 
of measuring the channel phase response to the BS. In 
this paper, we present an alternative method of achieving 
coherent transmission iteratively using a simple feedback 
control algorithm, which removes the need for explicit 
estimation of the channel to the BS, and greatly reduces 
the level of coordination required among the sensors. 

Other related work on synchronization in sensor net- 
works is based on pulse-coupled oscillator networks 
[8] and biologically inspired (firefly synchronization) 
[9] methods. These methods are elegant, robust and 
suitable for distributed implementation, however they are 
limited by assumptions of zero propagation delay and the 
requirement of mesh-connectivity, and are not suitable 
for carrier phase synchronization. 

We consider the following model to illustrate our 
ideas. The protocol is initialized by each sensor trans- 
mitting a common message signal modulated by a car- 
rier with an arbitary phase offset. (This phase offset 
is a result of unknown timing synchronization errors, 
and is therefore unknown.) When the sensors' wireless 
channel is linear and time-invariant, the received signal 
is the message signal modulated by an effective carrier 
signal that is the phasor sum of the channel-attenuated 
carrier signals of the individual sensors. At periodic 
intervals, the BS broadcasts a single bit of feedback 
to the sensors indicating whether or not the received 
SNR level increased in the preceding interval. Each 
sensor introduces an independent random perturbation 
of their transmitted phase offset. When this results in 
increased total SNR compared to the previous time 
intervals (as indicated by feedback from the BS), the 
new phase offset is set equal to the perturbed phase by 
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each sensor; otherwise, the new phase offset is set equal 
to the phase prior to the perturbation. Each sensor then 
introduces a new random perturbation, and the process 
continues. We show that this procedure asymptotically 
converges to perfect coherence of the received signals, 
and provide a novel analysis that accurately predicts 
the rate of convergence. We verify the analytical model 
using Monte-Carlo simulations, and use it to optimize 
the convergence rate of the algorithm. 

The rest of this paper is organized as follows. Sec- 
tion |n| describes our communication model for the 
sensor network. A feedback control protocol for dis- 
tributed beamforming is described in Section IIII-AI and 
its asymptotic convergence is shown in Section IIII-BI 
Section HV1 describes an analytical model to characterize 
the convergence behavior of the protocol. Some ana- 
lytical and simulation results are presented in Section 
Ivl Section IV-AI presents an optimized version of the 
feedback control protocol. Sections IV-Bl and fV-Cl present 
some results on scalability, and the effect of time-varying 
channels respectively. Section IVII concludes the paper 
with a short discussion of open issues. 

II. Communication Model for a Sensor 
Network 

We consider a system of N sensors transmitting a 
common message signal m(t) to a receiver. Each sensor 
is power constrained to a maximum transmit power of 
P. The message m(t) could represent raw measurement 
data, or it could be a waveform encoded with digital 
data. We now list the assumptions in this model. 

1) The sensors communicate with the receiver over 
a narrowband wireless channel at some carrier 
frequency, f c . In particular, the message bandwidth 
B < W c , where B is the bandwidth of m(t) and 
W c is the coherence bandwidth of each sensor's 



channel. This means that each sensor has a flat- 
fading channel to the receiver. Therefore the sensor 
i's channel can be represented by a complex scalar 
gain hi. 

2) Each sensor has a local oscillator synchronized to 
the carrier frequency f c i.e. carrier drift is small. 
One way to ensure this is to use Phase-Locked 
Loops (PLLs) to synchronize to a reference tone 
transmitted by a designated master sensor as in 
[6]. In this paper, we use complex-baseband nota- 
tion for all the transmitted signals referred to the 
common carrier frequency / c . 

3) The local carrier of each sensor i has an unknown 
phase offset, 7, relative to the receiver's phase 
reference. Note that synchronization using PLLs 
still results in independent random phase offsets 
7, = (2ttJ c ti mod 27r), because of timing 
synchronization errors n that are fundamentally 
limited by propagation delay effects. 

4) The sensors' communication channel is time- 
slotted with slot length T. The sensors only trans- 
mit at the beginning of a slot. This requires some 
coarse timing synchronization: r, <C T where Ti 
is the timing error of sensor i. 

5) Timing errors among sensors are small compared 
to a symbol interval (a "symbol interval" T s is 
nominally defined as inverse bandwidth: T s = i). 
For a digitally modulated message signal m(t), 
this means that Inter Symbol Interference (ISI) can 
be neglected. 

6) The channels hi are assumed to exhibit slow- 
fading, i.e. the channel gains stay roughly constant 
for several time-slots. In other words T s <C T <C 
T c , where T c is the coherence time of the sensor 
channels. 
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Distributed transmission model: The communica- 
tion process begins with the receiver broadcasting a 
signal to the sensors to transmit their measured data. 
The sensors then transmit the message signal at the next 
time-slot. Specifically, each sensor transmits: Sj(i) = 
A- gimit — Tj), where r, is the timing error of sensor 
i, A oc VP is the amplitude of the transmission, and g. t 
is a complex amplification performed by sensor i. Our 
objective is to choose gi to achieve optimum received 
SNR, given transmit power constraint of P on each 
sensor. For simplicity, we write hi = aie^ i and gi = 
bie> Bi , where bi < 1 to satisfy the power constraint. Then 
the received signal is: 

N 

r(t) =£toi(t)e™ +n(t) (1) 

i=l 

N 

= A^ /','/,' •'" '";/ - Tj) + n(t) 

i=l 
TV 

= aAe^+^+^M* - tj) + n(i). (2) 

i=l 

In the frequency domain, this becomes: 

iv 

i=l 

TV 

« A- M(f) £ aAe^ M ^ + N(f), (3) 
i=l 

where rc.(t) is the additive noise at the receiver and iV(/) 

is its Fourier transform over the frequency range | / 1 < 

B 

2 ■ 

In Q, the phase term ji accounts for the phase offset 
in sensor i. In (|3jl, we set e~ J ^ Ti w 1 because /r^ < 
Sri e|-<1. Equation motivates a figure of merit 
for the beamforming gain: 

N 

G=\Y^he^ +9 ^\ (4) 

i=l 

which is proportional to the square-root of received SNR. 

Note that fej < 1, in order to satisfy the power 
constraint on sensor i. From the Cauchy-Schwartz In- 



equality, we can see that to maximize G, it is necessary 
that the received carrier phases $^ = 7, + 9{ + ipi, are 
all equal: 

JV 

G= \Y,aihj**\ 

8=1 

N 

»=i 

iV 

= E(v^)(V5?^*')l (5) 

i=l 

JV 

< G opt = (E ! witri equality if and only if $ l = and 6j = 1 

i=l 

(6) 

However sensor i is unable to estimate either ji or ipi 
because of the lack of a common carrier phase reference. 
In the rest of this paper, we propose and analyze a 
feedback control technique for sensor i to dynamically 
compute the optimal value of Oi so as to achieve the 
condition for equality in 

III. Feedback Control Protocol 




Fig. 1. Phase synchronization using receiver feedback 

Fig- HI illustrates the process of phase synchronization 
using feedback control. In this section, we describe the 
feedback control algorithm, and prove its asymptotic 
convergence. 

A. Description of Algorithm 

The protocol for distributed beamforming works as 
follows: each sensor starts with an arbitrary (unsyn- 
chronized) phase offset 7$. In each time-slot, the sensor 
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applies a random perturbation to 9i and observes the 
resulting received signal strength y[n] through feedback. 
The objective is to adjust its phase to maximize y[n] 
through coherent combining at the receiver. Each phase 
perturbation is a guess by each sensor about the cor- 
rect phase adjustment required to increase the overall 
received signal strength. If the received SNR is found to 
increase as a result of this perturbation, the sensor adds 
the appropriate phase offset, and repeats the process. 
This works like a distributed, randomized gradient search 
procedure, and eventually converges to the correct phase 
offsets for each sensor to achieve distributed beamform- 
ing. Fig. |3 shows the convergence to beamforming with 
N = 10 sensors. 

iterations 10 iterations 




Fig. 2. Convergence of feedback control algorithm 

Let n denote the time-slot index and y[n] the ampli- 
tude of the received signal in time-slot n. From (0, we 
have: y[n] oc |X)i a i e ' 7 I where <&i[n] is the received 
signal phase corresponding to sensor i. We set the pro- 
portionality constant to unity for simplicity of analysis. 
At each time instant n, let f?j[n] be the best known 
carrier phase at sensor i for maximum received SNR. 
Each sensor uses the distributed feedback algorithm to 
dynamically adjust 0i[n] to satisfy (|6} asymptotically. 



The algorithm works as follows. 

Initially the phases are set to zero: 9i[l] = 0. At each 
time-slot n, each sensor i applies a random phase pertur- 
bation Si [n] to 0i [n] for its transmission. As a result, the 
received phase is given by: $j[n] = r Yi+8i[n]+Si[n]+il>i- 
The BS measures y[n] and keeps a record of the highest 
observed signal strength 0^, ea t[ro] = maxj <n y[k] in 
all previous timeslots. At the end of each timeslot, the 
BS broadcasts a one-bit feedback message that indicates 
whether the received signal strength of the preceding 
timeslot was higher than the previous highest signal 
strength. Depending on the feedback message, each 
sensor i updates its phase according to: 

9i[n + l] = l (7) 
I 6, [n] otherwise. 

Simultaneously, The BS also updates its highest received 
signal strength: 

y best [n + 1] = nmx(y best [n], y[n]j (8) 

This has the effect of retaining the phase perturbations 
that increase SNR and discarding the unfavorable ones 
that do not increase SNR. The sensors and the BS repeat 
the same procedure in the next timeslot. 

The random perturbation 5i[n] is chosen indepen- 
dently across sensors from a probability distribution 
8i[n] ~ fs(Si), where the density function fs(Si) is a 
parameter of the protocol. In this paper, we consider 
primarily two simple distributions for fs(Si): (i) the 
two valued distribution where Si = ±5q with proba- 
bility 0.5, and (ii) the uniform distribution where Si ~ 
uniform [— 60, So]. We allow for the possibility that the 
distribution fs(Si) dynamically changes in time. 

It follows from that if the algorithm were to 
be terminated at timeslot n, the best achievable signal 
strength using the feedback information received so 
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far, is equal to D4estM> which correspond to sensor i 
transmitting with the phase 9i[n]. 

3W[n] = Ea<e*'W| 

i 

where $ 4 [n] = 7, + S [n] + Vi (9) 

B. Asymptotic Coherence 

We now show that the feedback control protocol 
outlined in Section UlI-AI asymptotically achieves phase 
coherence for any initial values of the phases Let 3> 
denote the vector of the received phase angles We 
define the function Mag($) to be the received signal 
strength corresponding to received phase $: 

Mag(#) = 153^**1 ( 10 > 

i 

Phase coherence means $i = <t> 3 = Q CO nst, where 
^ const is an arbitrary phase constant. In order to remove 
this ambiguity, it is convenient to work with the rotated 
phase values fa — $i — $0, where $0 is a constant 
chosen such that the phase of the total received signal 
is zero. This is just a convenient shift of the receiver's 
phase reference and as dl 01 shows, such a shift has no 
impact on the received signal strength: 

Mag(^) = Mag($) (11) 

We interpret the feedback control algorithm as a 
discrete-time random process 34est W i n the state-space 
of fa the state-space being the TV-dimensional space 
of the phases fa constrained by the condition that the 
phase of the received signal is zero. We observe that the 
sequence {ybest [n] } is monotonically non-decreasing, 
and is upperbounded by G opt as shown in (|5}- Therefore 
each realization of {ybest [n] } is always guaranteed to 
converge to some limit Go < G opt . Furthermore Go = 
lim^oo y b est[n] = sup JWOO y bes t [n] i.e. the limit G 
of the sequence {ybest [n] } is the same as its Least Upper 



Bound [10]. In general the limit Go would depend on 
the starting phase angles fa We now provide an argument 
that shows (under mild conditions on the probability den- 
sity function fs(Si)), that in fact {34estM} converges to 
the constant G opt with probability 1 for arbitrary starting 
phases fa The following proposition will be needed to 
establish the convergence. 

Proposition 1: Consider a distribution fs(Si) that has 
non-zero support in an interval (— Sq,So), Given any 
^ ^ 0, and Ma,g(fa < G opt - e, where e > is 
arbitrary, there exist constants e x > and p > such 
that Prob(M&g((f> + 5)- Mag(^) > e x ) > p. 

Proof. For the class of distributions fs(Si) that we 
consider, the probability of choosing Si in any finite 
interval / C (—So, Sq) is non-zero. One example of such 
a class of distributions is fs(Si) ~ uniform[— <5 , <5 ]- 

Recall that the phase reference is chosen such that the 
total received signal J^. a i e^ i has zero phase. First we 
sort all the phases fa in the vector <f> in the descending 
order of | fa \ to get the sorted phases fa* satisfying | cj>\ \ > 
1 02 1 > ... > \ fa N \, and the corresponding sorted channel 
gains a*. We use the condition Mag(0) < G opt — e to 
get: 

«»(#) E < < E < cos (^) ^ G ° pt ~ e 

i i 

fa x >fa = ^- x { _ . ) (12) 

Now we choose a phase perturbation S\ that decreases 
\fa{\. This makes the most mis-aligned phase in 4> closer 
to the received signal phase, and thus increases the 
magnitude of the received signal. If fa > 0, then we 
need to choose a Si < 0, whereas if fa < 0, we need 
Si > 0. In the following, we assume that fa > and 
fa- > So. The argument below does not depend on these 
assumptions, and can be easily modified for the other 
cases. Consider Si G (— Sq,— This is an interval in 
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which fs(Si) is non-zero , therefore there is a non-zero 
probability p\ > of choosing such a <5i. We have: 



a* cos((?5* + Si) — a* cos(^>*) > 2ei 



where ei 



aT sin 



(13) 



We observe that ei and pi do not dependent on 0. 

The perturbation 8± by itself will achieve a non-zero 
increase in total received signal, provided that the other 
phases <f>* do not get too mis-aligned by their respective 

Si: 



Theorem 1: For the class of distributions fs(Si) 
considered in Proposition 1, starting from an arbitrary <f>, 
the feedback algorithm converges to perfect coherence of 
the received signals almost surely, i.e. y bes t [n] — ► G opt 
or equivalently cf>[n] — * (i.e. <pi[n] — > 0,Vi) with 
probability 1. 

Proof: We wish to show that the sequence ybest [n] = 
Mag(0[n]) — * G opt given an arbitrary <?5[1] = ^. 
Consider an arbitrarily small e > and define T e ((f>) 
as the first timeslot when the received signal exceeds 



G opt - e. 

Mag(<^ + S)- Mag(^) = ^ a*(cos(0* + Si) - cos^)) 

By definition if n < T £ (0), then J^estlT 1 ] 



Si) 



> 2e x + a? (cos( 



i>l 



- <5i) — co 
(14) 



^ag(^[n]) < G opt — e, and by Proposition 1, 

a*(cos((f)* + 8. t ) - cos(<f>*)) 
^roft(y fce;5t [n + 1] - ybestXA > £ i) > P f °r some 

nts ei > and p > 0. We have: 

E(y best [n + l]-y best [n]) >e lP yn<TS) (17) 



We note that since Mag(</>) is continuous in each of 
the phases <fi*, we can always find a > to satisfy: 



< (cos(</)* +Si)- C0S((7J*)) 



< 



Ar- 1 



,V|5,| < e, (15) 



In particular the choice = a ,^_^ , satisfies dl5> . and 
this choice of e t is independent of 0. With the <Vs chosen 
to satisfy d!5l >. we have: 



- £ i <E<( 

j>i 



cos 



+ -cos(0*)) < 



(16) 



Since fs{Si) has non-zero support in each of the non- 
zero intervals (— e;,£i), the probability p^ of choosing 
<5; to satisfy dl5t is non-zero, i.e. p^ > 0, which is 
independent of (f>. Finally, we recall that each of the Si 



Using ( TTtI we have: 

G opt >y best [n + l] 

n 

= y b est[l}+J2(ybest[k+l}-y b , 
k=\ 

n 

>^2(y b est[k + i}-y b e S t[k}) 

fe=l 

Taking expectation we have: 

n 

G o P t > E (j2(ybest[k + 1] - y best [k})) 

n 

?(X)(3W[* 



(18) 



> n -E 



> Prob(TS) 

> Prob[r t {(j)) > ;/ )//< ip 



l]-%e.t[* 
(19) 



are chosen independently, and therefore with probability where we obtained dl9l l by using dl7> . Therefore we 
p = J| i pi > 0, it is possible to find 8\ to satisfy dl 31 have Prob(T e ((j)) > n) < - — > 0,as n — > oo. 



and <5,j,i > 1 to satisfy dl5> . For 5 chosen as above 

Mag(<75 + 8) — Mag(0) > ei, and therefore Proposition shown that 3>f, es t[ra] 
1 follows. □ surely. 



Since this is true for an arbitrarily small e, we have 
G opt and 4>[n] -> almost 

□ 



T&) > n 
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IV. Analytical Model for Convergence 

The analysis in Section UlI-BI shows that the feedback 
control algorithm of Section IIII-AI asymptotically con- 
verges for a large class of distributions fs(Si); however 
it provides no insight into the rates of convergence. 
We now derive an analytical model based on simple, 
intuitive ideas that predicts the convergence behavior 
of the protocol accurately. We then use this analytical 
model, to optimize fs(di) for fast convergence. 




50 100 150 200 

timeslots, n 



Fig. 3. Motivating the Analytical Model: two simulated instances 
with N = 100, f s (Si) ~ uniform[- |j, 



A. Derivation of Analytical Model 

The basic idea behind our analytical model is that 
the convergence rate of typical realizations of y bes t [n] 
is well-modeled by computing the expected increase in 
signal strength at each time-interval given a distribution 
fs(Si). This is illustrated in Fig. [3] where we show 
two separate realizations of ybest [n] from a Monte-Carlo 
simulation of the feedback algorithm. 

We define the averaged sequence y[n] recursively as 
the conditional value of ybest [n + 1] given 34est [n] : 

y[l] = £(Mag(^[l])) (20) 

y[n + 1] = E m (y best [n + 1] y best [n] = y[n]) (21) 



The initial value y[X] in J20t is set under the assumption 
that the received phases cf>[l] are randomly distributed in 
[0,27r). For subsequent timeslots, ybest[n + 1] in (12 li 
is conditioned on 34estM but the phase vector <f>[n] is 
not known. Some remarks are in order regarding this 
definition, particularly the relationship of y[n] with the 
(unconditionally) averaged ybest [n] . Let 

y[n+l] = F(y[n]), where F(y) = E {y bes t[n+l]\y best [n] = y) 

(22) 

Consider: 



E{y b est [n + 1]) = Ey best[n] [E(y b6st [n + l}\y best [n}) 
= Ey bcst[n] (F(y b est[n})) 
~ F^E(y b est[n])) =y[n + l] (23) 

In most cases, the function F(y) is concave, and there- 
fore (by Jensen's Inequality) the approximation in d23l > 
represents an overestimate of the unconditional aver- 
age of ybest[n + 1]. Also in different instances of the 
algorithm, we would expect to see different random 
evolutions of cf>[n] and J-^est^] with time, and an av- 
eraged quantity only provides partial information about 
the convergence rate. Fortunately, as Fig. [5] shows, even 
over multiple instances of the algorithm, the convergence 
rate remains highly predictable, and the average charac- 
terizes the actual convergence reasonably well. Since the 
variation of the random ybest M around its average value 
is small, the approximation in d23l also works well. Our 
goal is to compute F(y) as defined in i22\ . 

Note that while i22\ is conditioned on ybest [n] being 
known, the phase vector <f>[n] is unknown. As 34estH 
increases, the phases fa [n] become increasingly clustered 
together, however their exact values are determined by 
their initial values, and the random perturbations from 
previous time-slots. In order to compute the expectation 
in J22L we need some information about 4>[n]. 
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We show in Section lTV-Bl that the phases 4>i[n] can be 
accurately modeled as clustered together according to a 
statistical distribution that is determined parametric ally 
as a function of J4est [n] alone. This is analogous to the 
technique in equilibrium statistical mechanics, where the 
individual positions and velocities of particles in an en- 
semble is unknown, but accurate macroscopic results are 
obtained by modeling the kinetic energies as following 
the Boltzmann distribution, which is fully determined 
by a single parameter (the average kinetic energy or the 
temperature). In our case, </>i[n] are modeled as indepen- 
dent and identically distributed (for all i) according to a 
distribution satisfying the constraint: 



] = Oi co s & = NE{a,)E 4>i (cos <j) t ) (24) 



y[n 



Therefore, even though the individual cj>i are unknown, 
we can compute all aggregate functions of <f> using this 
distribution, as if the fa are known. This is an extremely 
powerful tool, and we now use it to compute F(y) 
treating 4>[n] as a given. Section IIV-BI completes the 
computation by deriving the distribution used to specify 
4>[n] given y best [n]. 

From the condition 34estM = y[n], we have: 



y[n 



cos< 



(25) 



where we used the fact that the imaginary part of the 
received signal is zero for our choice of phase reference. 
We have the following expressions (omitting the time- 
index on <j)[n] and 8[n] for convenience): 



ybest [n + 1] = 



Mag(0 + <5) if Mag(<^) + S) > y[n] 



otherwise. 



(26) 



We now express Mag(0 + 6) as a sum of i.i.d. terms 
from each sensor, and invoke the Central Limit Theorem 



(CLT). 



(27) 



= | di (cos <pi cos d{ — sin (f>i sin 8{j + j a>i (cos (f>i sin S 

i i 

= | (C s y[n] +xi) + jx 2 , (28) 
where Cs = Eg (cos Si) , (29) 

di ^cos cj>i (cos Si — Cs) — sin <f>i sin SA , 

i 

(30) 

di ^cos sin Si + sin (pi cos 5A 

i 

(31) 

The random variables x\, x 2 are illustrated in Fig. |4] 

' Xd 



X\ 



x 2 




< 
< 



y[n]E(cosS) j Xl 

y[n] 



Fig. 4. Perturbation in the total received signal. 

Both X\ and x 2 as defined in J28t are linear com- 
binations of iid random variables, sin Si and cos Si. 
Therefore as the number of sensors N increases, these 
random variables can be well-modeled as Gaussian, as 
per the CLT [11]. Futhermore, Xi,x 2 are zero-mean 
random variables, and their respective variances a\,a\ 
are related by: 

ff i = \ E a i (C 1 - °?) ~ c°s(2^)(Q 2 - C 2S ) 

i 



where C 2 s = E s (cos(25i)) 



(32) 



With these simplifications, the statistics of y[n + 1] 
only depends on the density function fs(Si) through Cs 
and C 2 g. We have the following proposition. 
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Proposition 2: Assuming that the CLT applies for 
random variable X\, the expected value of the received 
signal strength is given by: 

y[n + 1] « y[n] (l -p-(l- C s )j + -±=e ^ 

(33) 

wh e re p = Q (&^.) (34) 

Proof. First we observe that the small imaginary 
component x 2 of the perturbation mostly rotates the 
received signal, with most of the increase in y[n + 1] 
coming from x\ (see Fig. |4}. 

Mag(0 + (5) = \C s y[n]+x 1 +jx 2 \ 

«(Cfa/[n]+xi) (35) 

Defining p as the probability that ybest[n + 1] > 
(I33> . (I34> readily follow from (I35> . (I26> using Gaussian 
statistics. □ 

We can rewrite d33l as: 

y[n + 1] = F{y[n]) = y[n] + f(y[n]) 

where f {y )±a ig (^^) 
1 

and = e 2 — xQ(x) (36) 

V 2vr 

Proposition 2 does not yet allow us to compute the 
y[n] because it involves the variance o\ that depends 
on the phases fa of the individual sensors. In the next 
section we present a statistical distribution for fa that 
allows us to calculate aggregate quantities such as o\ 
without knowledge of the individual fa . 

B. Statistical Characterization of Sensor Channels 

The statistical model is based on the assumption that 
each sensor has a channel to the BS of similar quality, 
and unknown phase. This means that the a/s are all 
approximately equal, and that the initial values of the 



phases fa are distributed independently 1 and uniformly 
in [0, 27r). In particular, we set a% = 1 for all sensors, 
which gives G opt = N. As the algorithm progresses 
towards convergence, the values of fa are distributed 
over a smaller and smaller range. In general, we expect 
that the distribution f<f,(fa) of fa[n] depends on the 
number of sensors N, the iteration index n, and the 
distribution of the perturbations f${8i). In the spirit of 
the statistical model, we consider large N, and look for 
a class of distributions that approximate f$(fa). 



0.25 




-2 -1.5 -1 -0.5 0.5 1 1.5 2 

Phase angle (radians) 

Fig. 5. Comparing a Laplacian Distribution with a Histogram of 
Empirically Observed Phase Angles 

We find that the Laplacian probability distribution 
gives the best results 2 in terms of accurately predicting 
the convergence behavior of the algorithm of Section 
IIII-AI Fig. [5] shows an empirically derived histogram 

'It is important to note that the fa are not random variables, however 
we statistically parametrize them using a probability distribution for the 
sake of compactness. 

2 The Laplacian distribution for fa is empirically found to work well, 
when compared with other families of distributions like the uniform 
and triangular distributions. 
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from a Monte-Carlo simulation of the feedback control 
algorithm. A Laplacian approximation is also plotted 
alongside the histogram. We now explain the details of 
the approximation. 

The Laplacian density function is given by [11]: 



1 -i^ii 
2<p 

For fa distributed according to d37l . we also have: 

1 



E [cos d>i) = 



E{cos24>i 



1 + 08 
i 



i 



(37) 

(38) 
(39) 



Therefore given that at iteration n of the feedback 
algorithm, the phase angles are <f>[n] = [<f>i<l>2»-<f>N]t we 
have: 



V 



[n] = Mag(0[n]) =^cne*+* = ^cos< 



(40) 



where we used a* = 1 in dlOl i. Now if we parametrize 
all the <pi using a Laplacian distribution, we can set <f>o 
such that Y^i cos 4>i = NEcj, (cos <f) . Thus we use J38i 
to rewrite ( I40l i as: 



y[n\ 



(41) 



We are now able to determine o\ given y[n]. 

Proposition 3: The variance a\ of X\ is given by: 



iV 



4-3 



JV 



(C£-tf M ) (42) 



Proof. Equation (1421 follows using d32l . and the 
value of the Laplacian parameter from (1411 along with 
the observation that ^Vcos(20j) = N ' E t j } (cos{2(j))) = 

□ 



N 



Using Propositions 2 and 3, we are able to analytically 
derive the average convergence behavior of the feedback 
control algorithm. In particular, we recursively calculate 
y[n] by substituting the variance a\ from d42l into d33l . 



C. Summary of Analytical Model 

We now summarize the analytical model derived in 
Sections IIV-AI and I1V-BI Our objective is to model 
the increase over time of the received signal strength 
by averaging over all possible values of the random 
perturbations. As mentioned before, we set the channel 
attenuations for each sensor to unity i.e. a,i = 1. 

1) Initially we set the received signal strength as 
y[l] = >/jV. This is the expected value of the 
signal strength if the initial phase angles are all 
chosen independently in [0, 2tt). 

2) At each time-interval (iteration) n > 1, given the 
probability distribution of the perturbations f$(8i) 
and the value of y[n], we compute the Laplacian 
parameter <fio using (I41t . and then compute the 
Gaussian variance a\ using (I42> and finally y[n + 
1] using the Gaussian statistics in d35l and i26\ . 

V. Performance Analysis of Feedback 
Control Protocol 

We now present some results obtained from the ana- 
lytical model of Section ffvl Fig. |6] shows the evolution 
of y[n) derived from the analytical model and also 
from a Monte-Carlo simulation with N = 100, for two 
different choices of the distribution fs(Si): a uniform 
distribution in [— gjj, gg] and a distribution choosing 
±-gjj with equal probability. The close match observed 
between the analytical model and the simulation data 
provides validation for the analytical model. 

We observe from Fig. [6] that the received signal grows 
rapidly in the beginning, but after y[n] gets to within 
about 25% of G opt , the rate of convergence becomes 
slower. Also while the simple two-valued probability 
distribution appears to give good results, it does not 
satisfy the condition for asymptotic coherence derived 
in Section IlILBl 



12 





simulated 
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500 1000 1500 2000 
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(a) fs(Si) ~ uniform(-~, ~) 



(b)/ 6 (di)~±^ 



Fig. 6. Comparison of Analytical Model with Monte-Carlo Simulation of Feedback Control Algorithm 



A. Optimizing the Random Perturbations 

In Fig. |6j we used the same distribution for the 
perturbations for all iterations of the algorithm. However 
this choice is not optimal: intuition suggests that it is best 
to choose larger perturbations initially to speed up the 
convergence and make the distribution narrower when 
the phase angles are closer to coherence. We now use the 
analytical model to dynamically choose the distribution 
fs(Si) as a function of y[n]. The general problem of 
choosing a distribution is a problem in calculus of vari- 
ations. Fortunately, it is possible to restrict ourselves to a 
family of distributions without losing optimality, because 
the analytical model only depends on the distribution 
through the two parameters Cg,C2S- Furthermore the 
parameters C§, C^g are highly correlated. To see this 
recall from < I3 1 1 and \32\ the definitions of Cg and C25 as 
the expected values of cos(£j) and cos(2<5i) respectively. 
Using the identity cos(2<5) = 2cos 2 ((5) — 1 and Jensen's 
Inequality, we can show that C28 is constrained by the 
value of Cg: 



2Ci - l<C 2 s<2C s -l 



(43) 



We are interested in Si corresponding to small random 
perturbations i.e. Si <§C 5. For such small values of Si, 
i43t allows only a small range of possible values of C25. 
Indeed we observe that cos(<5i) and cos(2Si) are very 
well approximated by the first two terms of the Taylor 
series: 

cos(*)»l-y, if (44) 

Equation <!44i indicates that both Cg and C25 are essen- 
tially determined by the second moment of Si, and there- 
fore even a one-parameter family of distributions fs(Si) 
is sufficient to achieve optimality of the convergence 
rate. Fig. [7] shows plots of the optimal choices of the 
(Cg, C 2 s) pair with N = 2000 over 10000 timeslots for 
two families of distributions: (i) the 3-point distributions 
P(±<5o) = p, P(0) = I — 2p parameterised by the 
pair (So,p), and (ii) the distributions uniformf— Sq, Sq] 
parametrised by Sq- At each iteration of the protocol, 
we used the analytical model to compute the value of 
the parameters (i.e. the pair (60, p) in case (i) and So 
in case (ii)) that maximizes the y[n + 1] given y[n]; 
the optimal parameters in each case were determined 
numerically using a simple search procedure. The two 
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curves in Fig.Qwere obtained by plotting (Cg, 6*2,5) pair 
corresponding to the optimal parameters for cases (i) and 
(ii) at each timeslot. The 3-point distribution is flexible 
enough to permit any (Cg,C2s) in the feasible region 
of J43i . For the example of Fig. it is clear that the 
uniform distribution achieves values of (C$, C25) that is 
close to optimal, thereby confirming the intuition of J44i . 
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Fig. 7. Near-Optimality of a One-Parameter Distribution 

We now use the family of distributions fs{Si) ~ 
uniform [—So, So] to obtain insight into the optimal con- 
vergence rate. Fig. [8] shows y[n] as a function of n for 
fixed values of So as well as for the optimized algo- 
rithm. We observe that the convergence rate decreases 
with time in all cases, and the optimized algorithm 
converges significantly faster than any fixed instance. 
Fig. [8] also shows the variation of optimal Sq with time. 
This confirms our intuition that at the initial stages of 
the algorithm, it is preferable to use larger perturbations 
(corresponding to large So), and when y[n] gets closer 
to G opt , it is optimum to use narrower distributions 
(corresponding to smaller So). 



Optimized algorithm compared to fixed fs(Si) 
iform[— <5n, So] for different <5n and N = 200 



B, Scalability Results 

We now turn to the analytical model to study the 
scalability of the feedback algorithm with the number 
of transmitting sensors N. We show the following scal- 
ability results: 

- The expected received signal strength at any time, 
always increases when more transmitters are added. 

- The number of timeslots required for the expected 
signal strength to reach within a certain fraction 
of convergence always increases with more trans- 
mitters, but increases no faster than linearly in the 
number of transmitters. 

Theorem 2: Let y\[n] and 2/2 M be the expected 
received signal magnitude at timeslot n when the number 
of transmitting sensors is N\ and iV2 respectively. If the 
sensors use the same distributions fs(Si) for all timeslots 
n, and N2 > Ni, then the following holds for all n: 



V2 [n] > yi[n] 



N 2 



(45) 
(46) 
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Proof. We offer a proof by induction. From Section 
HV^Cl we know that y 2 [l] > y x [l] and ^ > ^1. To 
prove J45i . we need to show that y 2 [n + 1] > yi[n + 1] 
given y 2 [n] > yi[n}. 

We write y x [n+ 1] = Fx (yx [n] ) , y 2 [n+ 1] = F 2 (y 2 [n] ) 
where -Fi(y) and are defined as in ( 1361 . Note that 

Fi(yi[n}) (and -FbCj/aM)) depends on the time index 
n not only through yi[n] (y 2 [n]), but also through the 
distribution fs(6i). We have suppressed this additional 
time-dependence to keep the notation simple. The func- 
tions F\(y) and i*2(y) satisfy the following properties: 

F 2 {y) >F 1 (y),Vy (47) 

Fi(y + ) > F 1 (y-), and F 2 (y+) > F 2 (y-) if y+ > y~ 

(48) 

To see this we observe from (1421 that for the same 
value of y, ci is larger for larger N, and since f{y) 
in (I36> increases with ax, ( I47> follows. To show (1481 . it 
is sufficient to show that F\(y) and ^2(2/) have a positive 
derivative with respect to y. This can be shown readily 
by differentiating the expression in d36l : 
dFx(y) d , _,_/y(l-C,)^ 



%(v+M) = i-(i-c s )Q( y[1 ai Us) ) >c s >o 



distributions. We apply Theorem 2 to the case where 
we use the distribution fs(Si) optimized for N x 
sensors in both cases. By definition y 2 [n] > y 2 [n], and 
yx[n] = yi[n], therefore y 2 [n] > yi[n],\/n. This proves 
J45t . Using the same argument for the distribution 
fs(Si) optimized for N 2 sensors, we can prove J46b . □ 

Another important criterion for scalability is the num- 
ber of timeslots Tf(N) required for the algorithm to 
converge to a fixed fraction, say / = 0.75 or 75% 
of the maximum for N transmitting sensors. Theorem 
2 shows that Tf(N) is an increasing function of N. 
Next we show that when the feedback algorithm is 
appropriately optimized, Tf(N) increases with N no 
faster than linearly. 

Theorem 3: The number of timeslots to convergence 
satisfies the following: 
Tf(N) 

lim — < tt, where tt is some constant. (51) 

Proof. First we use J43I to get a lower-bound for the 
variance a\. With y[n] = /■ G opt = /■ N we have: 



dy dy 

(49) 

We are now ready to complete the proof of ( 145 \ by 
induction. Given that y 2 [n] > y\[n], we have: 

y 2 [n+l] = F 2 (y 2 [n]) > F x (y 2 [ri\) > Fx(y x [n}) = yx[n-\ 

(50) 

where we used ( 147 \ and d48i for the two inequalities. 

This completes the proof of ( 145 \ . The proof of J46i 
by induction is similar and is omitted. □ 

Corollary: The scalability relations ( 145 \ and J46l > hold 

when the sensors use optimized distributions f${5i) in 
both cases. 

Proof. Let yi[n) and y 2 [n] be the expected received 
signal magnitudes using the respective optimized 



(1 - C s ) 2 < Cl - C 2S < (1 - Cf) 



(52) 



Using the upper bound from ( I52> in (I42> . we have: 



1 2 V4iV-3y[7i] 

1-/ 



> 2JV(1 - C s 



-4-3/ 

We now use a bound for the Gaussian Q-Function: 



(53) 



Q{x) > 



1 



« 2 A 



1 



/2tt x x 6 

Using ( I54> . we rewrite ( I33l > to get: 



(54) 



Ay[n] = y[n + 1] - y[n] > -^=e "=> (\ - \) 

V Z7T X X 



where x = 



y[n](i-c s ) 



(55) 
(56) 



15 



The bound in d55l has a maximum at xq w 3.6; 
choosing a C$ such that x is close to Xq, does not neces- 
sarily optimize the RHS in (I55> . because <j\ also depends 
on C$. However such a choice for Cg does provide a 
meaningful lower bound on the optimal Ay*[n]. 



2x Q ( l-f 
. ». . 2 / l-/w 1 -g. 1 



X Xq 



(57) 



(58) 



where I57i is obtained by backsubstituting d56l into ( I53l l. 
Let us denote the RHS of (58) by if (/). 

We observe that the lower bound in 1581 only depends 



on the fraction / = ^H. Let T fAf (N) be the number 
of timeslots required for the feedback algorithm to 
increase y[n] from a fraction / — A/ to a fraction / 
of convergence. If A/ is small enough, we can use d58i 
to write: 



C. Tracking Time-varying channels 

So far we have focused on the simple case of time- 
invariant wireless channels from each sensor to the BS. 
In practice, the channel phase response varies because of 
Doppler effects arising from the motion of the sensors 
or scattering elements relative to the BS. In the dis- 
tributed beamforming scenario, Doppler effects also arise 
because of drifts in carrier frequency between the local 
oscillators of multiple sensors. Therefore an important 
performance metric for the feedback control algorithm 
is its ability to track time-varying channels. Intuitively 
we expect that the algorithm should track well as long 
as the time-scale of the channel variations is smaller 
than the convergence time of the algorithm. In light 
of the scalability results in Section IV-BI the algorithm 
performs better for smaller N because the corresponding 
convergence time is smaller. 



Af-N = y[n]-y[n-T fAf (N)] 

« Ay[n]-T fAf (N) 
> K(f)T f , Af (N) (59) 
Therefore T fA f {N) < (60) 

Since T/ is just a sum of terms like 7/ ; a/> (15 1 i 
immediately follows. □ 



Theorem 3 is illustrated by the results in Fig.|9] where 
the number of timeslots required to get within a certain 
fraction of convergence is plotted against number of 
transmitters N for a fixed distribution (Fig. |9(a)} as 
well as optimized distributions (Fig. |9(b)) . These results 
show that the feedback algorithm is highly scalable with 
number of transmitters. 
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Fig. 10. Tracking Performance for Time Varying Channels: N = 100, 
Doppler rate=^y radians/timeslot 



A simulation of y[n) with time in the presence of 
channel time-variations is shown in Fig. [H)] This plot 
uses a fixed distribution for the phase perturbations, as 
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Fig. 9. Scalability of Feedback Control Algorithm with Number of Sensors 



the analytical model for optimization is not applicable 
to the time-varying case. A more detailed study of the 
tracking performance of the feedback control algorithm 
is beyond the scope of the present work. 

VI. Conclusion 

In this paper, we presented a simple algorithm for 
distributed beamforming in sensor networks, that is 
based on the idea of using SNR feedback from the 
receiver to perform phase synchronization in an iterative 
manner. This algorithm can be easily implemented in 
a decentralized manner and is guaranteed to achieve 
asymptotic coherence under mild assumptions. We also 
derived an analytical model that predicts the performance 
of the algorithm accurately, and offers insight into the 
convergence behavior. 

This paper represents an initial study into a new ap- 
proach to the problem of distributed synchronization, and 
leaves several open issues. We presented the Laplacian 
distribution to model the statistics of the phase angles 
<t>i as an empirical observation. However the underlying 
reason why the Laplacian distribution works so well 



is not clear. In addition the stability and convergence 
behavior of the feedback control algorithm under non- 
idealities like time-varying channels, and the effects of 
noise are open issues for future work. 

While we use the term "sensors" for the cooperating 
nodes performing distributed beamforming, the tech- 
nique developed here is of more general applicability. 
For example, it could be used as the basis for cooperative 
communication between clusters of nodes in a wireless 
ad hoc network. In such a context, it would be of interest 
to examine how the use of distributed beamforming 
would impact the design of medium access control and 
network layer protocols. 
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